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Based on a newly developed size-dependent cohesive energy formula for two-dimensional materials, a 
unified theoretical model was established to illustrate the gap openings in disordered graphene flakes, 
involving quantum dots, nanoribbons and nanoporous sheets. It tells us that the openings are essentially 
dominated by the variation in cohesive energy of C atoms, associated to the edge physicochemical nature 
regarding the coordination imperfection or the chemical bonding. In contrast to those ideal flakes, 
consequently, the gaps can be opened monotonously for disordered flakes on changing their sizes, affected 
by the dimension, geometric shape and the edge saturation. Using the density functional theory, 
accordingly, the electronic structures of disordered flakes differ to the ideal case because of the edge 
disorder. Our theoretical predictions have been validated by available experimental results, and provide us a 
distinct way for the quantitative modulation of bandgap in graphene for nanoelectronics. 

I n the fields of materials science, increasing interest is now being attracted to graphene, a two-dimensional (2D) 
I sheet composing of carbon atoms arranged in a honeycomb lattice structure 1 . Because of its compatibility with 
I industry- standard lithographic processing, electron mobility up to 150 times greater than Si and a thermal 
conductivity of graphene twice of that of diamond 2 , graphene shows its promise as a future material for the nano- 
and micro -electronics applications, such as field effect transistors (FETs) 3 ' 4 , biological membranes 5 , sensitive 
platforms 6,7 , and so on. Among them, as the most successful device concept in electronics, most works on 
graphene devices so far have been related to FETs 8 ' 9 . When graphene is used as the channel material, the 
graphene-based transistors will own an ultrafast speed, making it able to operate in the THz frequency range 10 . 
Moreover, the extremely thin channels will possibly allow the transistors to be scaled to shorter channel lengths 
without encountering the adverse short channel effects, which restrict the performance of existing devices. 

To achieve graphene-based electronic devices with high performance, specific bandgaps of graphene are 
usually required for different electronic applications 11 " 13 . As instance, since the I 0 JI 0 ff ratios in FETs are largely 
dependent on E g of the channel materials on the basis of their proportional relationship with I on /I 0 ff oc e Eg/ksT , 
graphene materials with sufficiently large bandgaps are demanded for room-temperature FET operation 14 " 16 , 
where E g denotes the bandgap, k B the Boltzmann constant and T the absolute temperature. In general, E g with 
only a few tens of meV is too weak for room-temperature operation of FETs 14 . To have I 0 JI 0 {{ ratios as 10-100, E g 
of graphene materials around 0.1 eV is required 15 , while for high I on /I 0 f{ of 10 7 , it should be large as 0.4 eV 1516 . 
Despite this, however, graphene cannot be directly adopted as semiconductors in electronics because of its zero 
bandgap [£ g (°°) = 0 eV] 17 , where oo denotes the bulk size in 2D. With such an issue, graphene creates a unique set 
of challenges for implementation into conventional electronics due to substantial leakage currents in the off state 
for field-effect devices. To address it, the degeneracy of the conduction and valance bands should be lifted via 
breaking the equivalence between the A and B sublattices 18 . Aiming at this purpose, the bandgap openings (BOs) 
in graphene have been realized by forming graphene quantum dots (GQDs) 2,19,20 , nanoribbons (GNRs) 2 ' 15 ' 21 " 23 and 
nanoporous graphene sheets (NPGs) 24 ' 25 . 

Upon the investigation on the BOs in GQDs, GNRs or NPGs, different edge natures should be considered. 
According to its honeycomb lattice, graphene may own ideal armchair (AC) or zigzag (ZZ) edges 12 . With the 
present fabrication techniques, however, these perfect edges are difficult to be formed for GQDs and GNRs 26 ' 27 . 
Instead, disordered edges are usually observed, which can be regarded as the mixture of AC and ZZ edge 
shapes 2 ' 28 . Because of it, the two edge structures of GNRs remain asymmetric 2 ' 26,27 . Moreover, their edges are 
kept naked with one dangling bond for each C atom. To stabilize those C atoms, dangling bonds can be chemically 
saturated with some radicals or groups denoted by R beneficial for the BOs, such as R = H 29 , F 30 or OH 31 . Such an 
edge saturation was not performed for GQDs and GNRs in most experiments due to its technical difficulty, while 
Ritter et at. 2 claimed that it can be performed through exposing GQD or GNR flakes to H radicals. Disordered 
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GQDs and GNRs without the edge saturation will be shown as GQD- 
null and GNR-null, while those saturated with R will be exhibited as 
GQD-R and GNR-R. As to mesoscopic NPGs 24 ' 25 , the holes in gra- 
phene are usually fabricated by the lithography technique with dif- 
ferent geometrical shapes, such as circles or squares (C-NPGs or S- 
NPGs). Far from the atomistic precision, the internal edges are usu- 
ally kept disordered and naked, referred to as DO-NPG-null. Owing 
to these attempts, the BOs were observed in disordered GQDs 19 ' 20 , 
GNRs 315 ' 21 ' 23 ' 32 " 34 and NPGs 24 ' 25 , while the D-dependence of E g (D) 
varies considerably among them, where D denotes the diameter of 
GQDs or holes in C-NPGs, width of GNRs, or the side length of holes 
in S-NPGs. On lowering their D, a monotone increase in the 
openings was observed for disordered and naked GQDs and 
GNRs 31519 " 21 ' 23 ' 32 " 34 , while it is more substantial for GQDs. As regards 
NPGs, a monotone decrease in the BOs has been reported 24 . In 
contrast, the BOs are lowered when their edges are saturated 2 . 
However, our understanding to the openings in disordered flakes 
is still in incubation. 

Much effort has been devoted to investigating the opening mech- 
anism in graphene. As a powerful technique, the computer simu- 
lation approaches, such as the first principles density functional 
theory (FP-DFT), Dirac equation and the tight-binding (TB) approx- 
imations, have been adopted for such studies. It may provide us band 
structure images, electronic density and other basic information of 
graphene for analyzing and explaining the electronic performance 35 . 
By this means, the opening behavior is usually studied using ideal 
GQDs and GNRs with regular edge structures 13 ' 29,31 ' 36,37 or ideal 
NPGs with perfect periodic antidot superlattices or holes 38 " 40 . With 
these contributions, however, it seems still difficult to illustrate the 
openings in disordered flakes, since these ideal flakes exhibit much 
different D-dependences of E g (D). The BOs in ideal GQDs and GNRs 
are relevant to edge structures. AC-GQDs have larger E g (D) than ZZ- 
GQDs do 36 . By exploring the electronic properties for a wide class of 
nanoribbons, GNRs exhibit a rich variety of bandgaps from metals to 
typical semiconductors 13 ' 37 . Typically, ZZ-GNRs show a metallic or 
semiconducting state 13 ' 29 ' 31 ' 37 , while AC-GNRs present the oscillating 
BO behaviors from the metallic to semiconducting state with respect 
to D 13,29 ' 37 , where the metallic points emerge at the period of 3 with iV a 
= 3p + 2, N a the number of AC chains and p a positive integer. 
Concerning ideal NPGs 38-40 , the opening behaviors are much com- 
plicate, which may be affected by not only the hole diameter but also 
the hole geometry, lattice geometry, separation between adjacent 
holes, hole locations and the edge geometry. When D of holes rises, 
E g (D) would increase for some antidot lattices 39 , such as hexagonal or 
triangular antidot lattices with a hexagon -shaped hole, while for the 
rectangular, rotated triangular and the honeycomb antidot lattices 
with a hexagon -shaped hole, only particular hole separations lead to 
an obvious opening 39 . The rectangular antidot superlattices perfo- 
rated with four holes of honeycomb symmetry can exhibits the semi- 
conducting or metallic states, depending on the width ratio of AC 
chains to ZZ chains 38 . More impressively, in particular, AC- 
GNRs with iV a = 1 1 can be changed from the metallic state to the 



semiconducting state when it is perforated with a honeycomb hole 
especially at the central position, while it is still metallic if a triangu- 
lar-shape hole is introduced 38 . 

Besides the above studies, some other theoretical ways have been 
explored for the BOs in DO -GNRs by considering the impact from 
the edge disorder 31 . With the renormalized lateral confinement 
approach 41 , E g (D) of DO-GNRs is decided by the structural disorder 
at the edges caused by etching or inaccuracies in the width deter- 
mination, while such a disorder can result in an "effective transport 
width" smaller than the ribbon's nominal width. Alternatively, the 
BO behaviors can be explained with the Coulomb blockade or char- 
ging effects 22 ' 26 ' 27 , where the BOs are induced not only by confine- 
ment but also by Coulomb blockade originated from the edge 
roughness. However, the BOs in disordered GQDs and NPGs have 
been less discussed. As shown in Table 1, in addition, some expres- 
sions have been presented to plot E g (D) of disordered GQDs, GNRs 
or NPGs. In general, these expressions are achieved using the fitting 
method with those experimentally measured E g (D) values, while one 
or two adjustable parameters a and/or b are present there. Since these 
parameters are case- dependent and vary considerably among differ- 
ent researches, the reported fitting expressions cannot take their 
essential roles in helping us revealing the opening mechanism either. 

Due to the limitation in existing approaches or expressions, some 
issues on the BOs in disordered graphene flakes have not been 
addressed. As instances, a theoretical way to describe the difference 
in E g (D) induced by the dimension effect is still absent. It also 
remains unclear how to depict the reverse change in E g (D) of 
NPGs relative to those of GQDs and GNRs. As to the role from 
the edge saturation, on the other hand, its exact mechanism has 
not been clarified as yet. On account of it, one still wonders how 
the BOs in saturated flakes will be weakened relative to those naked 
ones. In practice, the gap opening mechanism of disordered flakes is 
much complicated, which is relevant to not only the size but the 
dimension, geometric shape and the edge saturation. However, a 
systematic theoretical investigation on their roles is still absent, 
which obstructs our physical insight into the BOs and throws obsta- 
cles in our way toward the full-scale application of graphene in 
electronics. Because of these, a unified theoretical approach should 
be developed to elucidate the BOs in disordered graphene flakes. 

In recent years, with the use of the mesoscopic thermodynamics of 
materials 42 ' 43 , a bandgap thermodynamics has been proposed to illus- 
trate the E g (D) function of 3D III-V and II-VI nanosemiconductors 
with macroscopic parameters of materials based on the nearly- free- 
electron approach 44 46 . In light of it, strikingly, it reads, 

E g (D)=[2-E c (D) 3D /E c (co)]E g (co). (1.1) 

In equation 1.1, £ C (D) 3D is the energy required to break the atoms 
in the solid into isolated atomic species, which can be obtained 
experimentally by measuring the latent heat of sublimation at any 
temperature, and extrapolating to zero Kelvin. Inspired by 
Lindemann's criterion for solid melting and Mott's expression for 



Table 1 Equations or expressions reported to depict measured E g (D) of disordered GQDs, GNRs and NPGs 


with or without the edge 


saturation 










flakes 


saturation 


E g (D) in eV 


a 


b 


GQDs 


H 


E Q (D) = a/D b * 


1.57 ±0.21 


1.19 ± 0.15 


GNRs 


null 


E g (D) = a/D 15 ' 21 


0.8eVnm 15 










1.26 eVnm 21 








E g (D) = a/(D - k) 23 41 


0.38 eVnm 23 


1 6 nm 








0.2 eVnm 41 


1 6 nm 






E g (D) = a/(De bD ) 22 23 


1 eVnm 22 


0.023/nm 








2 eVnm 23 


0.026/nm 


NPGs 




E g (D) = o/(33 - D) 24 


0.95 eVnm 
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vibrational melting entropy 47 " 50 , utilizing the proportional relation- 
ship between E C (D) and T m (D), £ c (D) 3D AE c (o°) function for the 3D 
nanocrystals having free surfaces with broken bonds is given as 44,45 , 

E c (D) 3D /E c (co)= exp[-(a° D -l)/(D/D» D -l)], (1.2) 

where oc® D is a surface physicochemical amount with oc® D = 
a sy (D) 2 / a in _ 3D (D) 2 , a 2 the mean square displacement of thermal 
vibration at the melting temperature, and the subscripts sv and in- 
3D the surface and interior atoms. Using Lindemann's criterion for 
solid melting and Mott's expression 49,50 , it was derived that o^ D = 
2S vi b( 00 )/3 J R + 1, where S vib (o°) denotes the atomic vibration entropy 
and R the gas constant. D/D® D = v/sh or D® D = 2(3 — d)h denotes the 
critical diameter of crystals where almost all atoms or molecules are 
located on the surface and a crystalline structure is no longer stable, 
v/s means the volume/surface ratio, d the dimension of the crystal 
with d = 0 for nanoparticles, d = 1 for nanowires and d = 2 for thin 
films. In light of equations 1.1 and 1.2, the bandgap expansion for 
nanoscaled 3D III-V and II-VI semiconductors induced by the 
surface or interface effect has been illustrated successfully 44-46 . 
Stimulated by this success, one may conceive that the aforesaid issues 
on the BOs in the 2D graphene system may also be addressed in the 
thermodynamic way. Figure 1 illustrates the schematic plots of (A) 
GQDs, (B) GNRs, (C) C-NPGs and (E) S-NPGs, while a tripolar flake 
marked with YA^Bx B2QC2 in (D) and a cross flake with 
+ A 1 A 2 B 1 B 2 C 1 C2D 1 D2 in (F) denote the nanoscopic functional units 
of C-NPGs in (C) and S-NPGs in (E), respectively. In view of it, an 
atomic layer can be observed at the edge of GQDs, GNRs and tripolar 
or cross flakes, which should be responsible for the BOs in graphene. 
With this way, the distinct edge physicochemical nature of C atoms, 
such as coordination imperfection (CI) or chemical bonding, should 
play an important role in the BOs in graphene with relation to the 
dimension, geometric shape and the edge saturation. However, no 
attempts have been made on it so far. 





(A) GQDs 

^graehene ^ 



(B) GNRs 





(C) C-NPGs 



(D) tripolar unit of (C) 





(E) S-NPGs 



(F) cross unit of (E) 



Figure 1 | Illustration of edge atomic layers (marked in red) for (A) 
GQDs, (B) GNRs, (C) C-NPGs, (D) tripolar unit of (C), (E) S-NPGs, and 
(F) cross unit of (E). 



3D ^ 

(NJN in -3DJ 
and iVin-2D denote the CN in the 



In this paper, a general thermodynamic approach will be proposed 
to illustrate the BO behaviors in disordered graphene flakes. Based on 
this, significantly, it provides us a new physicochemical insight into 
the openings in graphene. With the DFT simulation, accordingly, the 
electronic structure of disordered GNRs will be investigated. 
Suggestions will be made how to tune E g (D) of graphene for its 
application in nanoelectronics. 

Results 

BOs in disordered and naked GQDs, GNRs and NPGs. In light of 
equations 2.1 and 2.2 in methods, E g (D) of naked graphene flakes can 
be worked out as long as oc 2D and D 2D are known. Since ol® d and oc 2D 
are decided by the CI of atoms at the surface of 3D nanocrystals or 
edge of 2D graphene flakes, oc 2D can be correlated with at® D via the 
coordination number CN of atoms at the surface N sv and that at the 
edge iV edge . According to the inverse proportion relationship between 
a 2 and CN of atoms, based on the definitions of oc° 3D in Introduction 
and %2D i n Methods, one has oc® D oc (iV sv AZV in _ 3D ) _1 and oc 2D oc 
(N edge /N in - 2D )-\ where N ir 

bulk 3D and 2D crystals. It thus reads oc^/a^ = (N sv IN in - 3I) )l 
(iV edge AZV in _ 2D ). Since iV in _ 3D = 12 and N sv = 9 for face-centered 
cubic crystals and iV in _ 3D = 8 and N sv = 6 for body-centered cubic 
ones 47 , N sv /N in - 3D = 3/4. Moreover, iV edge = 2 and iV in _ 2 D = 3 for 
2D graphene flakes. On the basis of these, oQ D can therefore be 
elucidated with oc 2D = 90^/8. Utilizing the published oc® D expres- 
sion cited in Introduction, one has = 9[2S vi b( 00 )/3 J R + l]/8. 

As regards D 2D of GQDs, GNRs and NPGs, it can be given with 
equation 2.4. Since lis = 2(2 - d)ID, D° 2D of GQDs and GNRs can be 
given as D/D° 2D = D/[2(2 - d)h] or D° 2D = 2(2 - d)h, where d = 0 for 
GQDs and d = 1 for GNRs. As regards NPGs, D 2D cannot be decided 
directly with that for GQDs or GNRs. This is because D here denotes 
not the size of functional tripolar or cross flake but that of holes in 
NPGs. Given that the effective size of tripolar or cross flakes is D eff , 
D e ff/D 2D should be explored using D e ff/D 2D = sllh. For the tripolar 
unit of C-NPGs in Figure 1(D), the edge-C atoms along the concave 
edges of A!A 2 , BiB 2 and QC2 own broken bonds, while those along 
straight lines of A 2 B l5 B2Q and C 2 Ai are still embedded in graphene 
without any CI, where not the latter but the former takes effect in the 
BOs of NPGs. Hence, the effective free edge / of YA 1 A 2 B 1 B 2 C 1 C 2 
comes from three concave sections of A^, BiB 2 and QC2 only. Let 
the distance between the center of any two neighboring circles O l5 0 2 
and 0 3 be C larger than D, I = A X A 2 + B^ + QC2 = 3A!A 2 = nDI 
2, and s = \/3C 2 /4 - 7tD 2 /8. Thus, we have lis = 4ttD/(2v / 3C 2 - 
tlD 2 ). Taking it into the above D e fflD 2D expression, one has D e fflD 2D 
= (2 V^C 2 - nD 2 )l(4nDh). Concerning the cross unit of S-NPGs in 
Figure 1(F), in analogy, lis = 4D/(C 2 - D 2 ), and D e{{ ID 0 2D = (C 2 - 
D 2 )l(4Dh). Inserting these D ef flD 2D expressions into equation 2.2 
instead of DID 2D , E g (D) as the function of D can be explored for 
C- and S-NPGs. 

With the above amounts, Figure 2 shows the solid E g (D) curves as 
the function of D for (A) GQD-null or GNR-null, (B) C-NPG-null 
and (C) S-NPG-null plotted with equations 2.1 and 2.2 . For the 
comparison purpose, the reported reciprocal or exponential fitting 
expressions are also shown using dashed and/or dotted curves for 
GNR-null 15 ' 41 ' 21 ' 23 and NPG-null 24 . Available experimental results 
(symbols) were shown accordingly for verification. In light of 
Figure 2, the bandgaps are substantially opened for naked GQDs, 
GNRs and NPGs. As D shrinks, E g (D) of GQDs and GNRs is 
expanded exponentially. Moreover, GQDs have larger E g (D) values 
than GNRs do. In contrast, much similar D- dependence of E g (D) has 
been observed for C-NPGs in (B) and S-NPGs in (C), where E g (D) 
shrinks on decreasing D, and E g (D) becomes strong as C declines 
from 44 nm to 22 nm. Our predicted curves go through experimental 
results shown in Figure 2 (A) -(B), supporting that our theoretical 
approach can be utilized to elucidate the BOs in graphene regarding 
the dimension and geometrical shape. As regards the upward or 
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downward deviation of reported E g (D) values of GNRs from our 
predictions, it should be attributed to the methods to measure 
E g (D). As exampled by GNRs, in particular, its E g (D) is usually 
derived from the I 0 JI o{{ ratios of GNR-based FETs 3 ' 15 ' 16 ' 51 . Utilizing 
this means, there will be some possible errors caused by inaccurate D 
determination due to over- etching underneath the etch mask, uncer- 
tainties in the D measurements based on atomic form microscope or 
the assumption of negligible Schottky barrier for holes in ultra-nar- 
row GNR-based FETs 15 ' 41 . 

BOs in disordered but saturated GQDs and GNRs. Using 
equations 2.1 and 2.2 , E g (D) of saturated GQDs and GNRs can be 
explored if is known, where D 2D is the same as that of naked 
GQDs and GNRs. Since oc 2 d is related to the edge atomic nature, a^, 
may be resolved via a^. According to equation 2.3, it reads oc 2D = 

^ed g e(^) 2 /^in-2D(^) 2 and a° 2D = a° edge (D) 2 / ' a in _ 2D (D) 2 . Divided 

*2D by a° D , it has a* D /a° D = ^ dge (D) 2 /^ dge (D) 2 . Since E c oc a" 2 , 

one gets a* D = [E° c (D)/E*(D)]ol° 2D , where E° C (D) and Ef(D) are the 
respective atomic cohesive energies of edge-C atoms of naked and 
saturated GQDs and GNRs. Provided that E° C (D) and Ef(D) have the 



(A) DO-GQD/GNR-null 

GQDS this work 

ref.15 

ref. 21 




Figure 2 | The solid E g (D) curves as the function of D using equations 2.1 
and 2.2 for (A) DO-GQD-null and DO-GNR-null, (B) DO-C-NPG-null 
and (C) DO-S-NPG-null at C = 22, 33 and 44 nm, where ot° 2D and D° 2D are 
from Results. For the comparison purpose, the reported reciprocal or 
exponential fitting expressions from Table 1 are also shown using dashed 
and/or dotted curves with E g (D) = a/D 15 > 2 \ E g (D) = a/(D - b) 23A1 and 
E g (D) = a/(De to ) 22 ' 23 forGNR-nullin (A) and£ g (D) = a/(33 - D) 24 in (B). 
The symbols denote available measured results with T 19 and ^ 20 for GQDs 
and D 41 , iV 5 , © 21 , ^ 23 , O 3 , A 32 , V 33 , and O 34 for GNRs in (A) and O 24 for 
C-NPGs in (B). For calculation, h = 0.142 nm 18 ' 35 . S vib (°o) is given with 
S v ib(°°) ~ S m (°°) — i? 44 ' 45 ' 50 , where S m (°°) means the melting entropy using 
S m (°°) = H m (oo)/T m (oo) 44 > 50 with H m (oo) and T m (oo) being the melting 
enthalpy and melting point, respectively. S vi b (°°) or S m of graphene can be 
referenced to that from graphite with H m (°°) = 132 kj/mol and 
T m (oo) = 5500 K 58,59 . 



same D- dependences, the ratio of E®(D)/E^(D) should be D- 
independent, and thus a^, = [£^/^] a 2D- With its ne lp> a m can ^ e 
given when E® and E^ are known. 

With the above derived expressions, Figure 3 shows the solid E g (D) 
curves as the function of D in light of equations 2.1 and 2.2 for DO- 
GQD-R and DO-GNR-R with (A) R = H, (B) R = F and (C) R = 
OH. In comparison, the reported reciprocal fitting expression 2 to 
depict the D-dependences of E g (D) for GQD-H as aforementioned 
is also given using dashed curves in (A). Available experimental 
results (symbols) were also shown accordingly. As can be seen in 
Figure 3, when the edge saturation process is carried out, the BOs are 
realized, where E g (D) is enhanced as D decreases. There exists a fair 
agreement between our curves and experiment results in Figure 3A, 
which validates our theoretical models. 

Investigation on electronic structure of DO-GNRs. Because of its 
importance in electronics, in usual, ideal AC or ZZ-GNRs have been 
adopted to investigate the electronic structure of graphene-based 
nanostructures 13,18 ' 29,37 ' 52 . However, seldom attempts have been 
made to it for DO-GNRs. Hence, the electronic structure of DO- 
GNRs with the AC or ZZ skeleton along the periodic axis, referred 
to as DO-GNRs(AC) or DO-GNRs(ZZ), is investigated here. Figure 4 
shows the band structures and the corresponding charge density 
isosurfaces of the lowest unoccupied molecular orbital (LUMO) 
and the highest occupied molecular orbital (HOMO) at T point of 




0.0 r ^-r~i--r-T^F=F=F^ 

0 5 10 15 

D (nm) 



Figure 3 | The solid E g (D) curves as the function of D in light of equations 
2.1 and 2.2 for DO-GQD-R and DO-GNR-R with R = H in (A), F in (B) 
and OH in (C), where a^, and D\ B are from Results. In comparison, the 
reported reciprocal fitting expression from Table 1 to depict the D- 
dependences of E g (D) with E g (D) = alD h for DO-GQD-H 2 is also given in 
(A) using dashed curve. The symbols denote the measured results with ■ 2 . 
For naked GQDs or GNRs, E° c = E c=c where E c=c is the bond energy of C 
= C. Upon the edge saturation with R, Ef = E c=c + £ c _ R /2. E c=c and 
£ C _ R will be taken as the averaged values from ZZ and AC-GNRs. The 
respective E c =c values are taken as 4.68 and 5.81 eV for £° of ZZ and AC- 
GNRs 56 . To have the respective E c=c and £ C -r are 4.97 eV and 5.47 eV 
with R = H, 4.60 eV and 5.22 eV with R = F and 4.44 eV and 5.06 eV with R 
= OH, where the values of these parameters are given with the simulation 
method in Supplementary Information. See the caption of Figure 2 for 
other necessary parameters in calculation. 
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DO-GNRs(AC) and DO-GNRs(ZZ) in (a) R = H, (b) R = F and (c) 
R = OH or of DO-GNRs(ZZ) in (d) R = H, where their width is 
given at D ~ 1.2 nm. In light of their band structures, the energy 
bands for spin-up and spin-down states are degenerate, while the 
openings are observed substantially for all DO-GNRs. The respective 
gaps are 0.44, 0.50 and 0.56 eV for DO-GNRs(AC) saturated by R = 
H, F and OH. Regarding DO-GNRs(ZZ) saturated with R = H, E g (D) 
= 0.40 eV, close to that of DO-GNRs(AC) with R =H. In 
comparison, some differences in HOMO and LUMO can be 
observed between DO-GNRs and ideal GNRs. In ideal AC-GNRs, 
the charge densities are concentrated homogenously along the C=C 
bonds over the whole framework 36 . In Figure 4(a)-(c) for DO- 
GNRs(AC), the HOMO and LUMO orbitals distribute mainly 
along the C=C bonds. However, the densities are inhomogeneous 
over the whole region, much weak in the inner area instead. 
Moreover, some electronic states appear at single ZZ-structured 
edge-C atoms (labeled as 1, 2 and 3). In ideal ZZ-GNRs, on the 
other hand, the charge densities are centered at single C atoms 
especially for those along the two ZZ edges 36 . As for DO- 
GNRs(ZZ) in (d), however, the LUMO and HOMO charge 
densities are mainly centered at single C atoms, which seems 
strong in the interior area. Some charge densities along the C=C 
bonds have also been observed. 



Discussion 

According to our predictions, the BOs observed for disordered and 
naked GQDs, GNRs or NPGs in Figure 2 are contributed from the 



decrease in E c of C atoms at the edge compared to that in the interior. 
In a perfect graphene, each C atom has a planar sp 2 configuration 
with three C = C bonds. When graphene flakes are naked, each edge- 
C atom has one dangling bond, and it thus consists of two C=C 
bonds only, leading to a decrease in E c of edge-C atoms relative to 
those in the interior. Consequently, the equivalence between the A 
and B sublattices is broken, giving rise to the BOs in graphene flakes. 
As regards the saturated case in Figure 3, the dangling bonds of edge- 
C atoms are terminated, while each edge-C atom owns two C = C 
bonds and one C-R bond, different from the case in the interior. On 
account of it, the equivalence of A and B atoms is broken, resulting in 
the BOs in GQD-R or GNR-R. This infers that the BOs in saturated 
GQDs and GNRs should also be ascribed to the bonding change of 
edge-C atoms. With regard to the difference in E g (D) among GQDs, 
GNRs and NPGs, it agrees with their D-dependences of lis ratio, 
which rises for GQDs and GNRs but shrinks for NPGs on lowering 
D (Supplementary Fig. SI). Concerning the large E g (D) of GQDs 
relative to GNRs, it should be attributed to the observation that 
GQDs have larger Us ratios than GNRs do. It should be noted that 
the D-dependences of E g (D) for GNR-null in Figure 2(A), NPG-null 
in Figure 2(B) and GNR-H in Figure 3(A) can be reflected on 
decreasing D with those reported fitting expressions. However, these 
expressions can only be adopted for elucidating the BO behaviors in 
those graphene flakes numerically, while neither of them can be 
adopted to address the issues on the essential roles played by the 
dimension effect, geometric shape and the edge saturation, ascribed 
to the presence of adjustable parameters a and/or b in them. In 
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Figure 4 | Band structure, charge density isosurfaces of LUMO and HOMO at T point of DO-GNRs(AC) with (a) R = H, (b) R = F and (c) R = OH and 
of DO-GNRs(ZZ) with (d) R = H. The E F is set to zero. Blue and yellow denote the positive and negative wave function contours, respectively, and the 
value of the isosurfaces is 0.02 e/A 3 . C, H, F and O atoms are denoted by gray, white, cyan and red balls, respectively. 
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Figure 2(A), in particular, several fitting expressions have been 
reported to elucidate the D- dependence of E g (D) for GNR-null. 
However, they vary considerably among different researches, and 
it is difficult for one to see which one is suitable for illustrating the 
D-dependence of E g (D) for naked GNRs. 

According to the previous studies, a monotonous increase of E g (D) 
has been reported for ideal AC and ZZ-GQDs 36 , which agrees with 
our predictions on E g (D) of DO-GQDs. Concerning the difference in 
the size of E g (D) between AC-GQDs and ZZ-GQDs 36 , it should be 
attributed to their different edge natures 53 . However, the opening 
behaviors of disordered GNRs differ to those ideal AC and ZZ- 
GNRs. The oscillating BO behaviors in AC-GNRs are induced by 
the full Fermi wavelength effect 13 . As for ZZ-GNRs, the reported 
metallicity is caused by the flat bands at the Fermi level E F from 
the highly localized states formed at the edges 18 , while the semicon- 
ducting state is stemmed from the staggered sublattice potential on 
the hexagonal lattice induced by the edge magnetization 29 . One may 
thus wonder whether the roles from the Fermi wavelength effect, the 
flat bands at E F or the interedge magnetic interaction observed in AC 
or ZZ-GNRs would affect the BOs in disordered GNRs or not. 
Among them, in fact, the roles from the Fermi wavelength effect 
and the flat bands at E F in DO-GNRs should be negligible. This is 
reasonable since almost no oscillating BO behaviors or metallic states 
have been reported for experimentally fabricated DO-GNRs, while 
the openings in them are observed substantially and rise monoto- 
nously on lowering D instead 315 ' 21,23 ' 32 " 34 . As to the role from the 
interedge magnetic interaction, it was investigated in the present 
work by performing spin-polarized and non spin-polarized calcula- 
tions for DO-GNR-R shown in Figure 4 (a) -(d). The obtained results 
show that almost no energy difference can be found for them 
between these two methods, indicating that all DO-GNRs give the 
nonmagnetic characteristics. In light of it, the contribution from the 
interedge magnetic interaction to the BOs in DO-GNRs should also 
be neglected. Obviously, the negligible roles from these factors 
should be attributed to the edge disorder. In Figure 4, in addition, 
the increasing order of E g (D) for DO-GNR-R at D ~ 1.2 nm from the 
simulation method with R = H, F and OH is the same as that from 
our predictions in Figure 3. However, their sizes are somewhat lower 
than those from our predictions, which might be related to the 
investigation approaches. As to the difference in the HOMO and 
LUMO charge density between DO-GNRs and ideal GNRs, it should 
also be relevant with the edge disorder in the former. 



As mentioned above, in addition, the BOs contributed from the 
confinement role of D of holes have been reported for some typical 
ideal antidot superlattices 38 " 40 , which are in agreement with our 
predictions on DO-NPGs. In contrast, the metallic state observed 
in some other superlattices should be induced by the flat bands 
attributed to some factors, such as the instability of the incomplete 
benzenoid bonding pattern, large width of the ZZ chains or the 
presence of the ZZ rims 38 " 40 . As the metal-semiconductor transition 
in porous AC-GNRs with N a = 1 1 regarding the position effect, it is 
induced by the suppression of the localized states 38 . However, it 
seems that the roles from these factors can be neglected in the BOs 
in DO-NPGs, since they can be much suppressed by the internal edge 
disorder of holes. As evidence, no metallicity has been observed in 
experimentally fabricated DO-NPGs 24 ' 25 ' 54 . Besides it, a more recent 
investigation on the BOs in ideal AC or ZZ-GNRs perforated with 
disordered holes claims that the sensitivity of E g (D) to the lattice 
symmetry is significantly suppressed in the presence of disorder 55 . 

To see how E g (D) of saturated DO-GQDs or GNRs will evolve 
relative to the naked case, Figure 5 shows AE g (D) as the function of D 
with AE g (D) = E g (D) R — E g (D) nu u. AE g (D) are negative for the edge 
saturation by H, F and OH, decreasing on the decline of D. This 
suggests that the BOs are weakened after the edge saturation. In 
comparison, the D-dependence of AE g (D) for GQDs is stronger than 
that for GNRs where the D-dependence of AE g (D) is the strongest for 
the H saturation and the weakest for the OH saturation. In view of it, 
further, one can see that E g (D) of GQDs or GNRs can be lined in 
sequence with E g (D) H < £ g (D) F < E g (D) OH < £ g (D) null . Frankly, 
such a variation is essentially relevant to the difference in E c of edge- 
C atoms for GQDs or GNRs with or without the edge saturation. 
Relative to naked graphene flakes, E c of edge-C atoms is increased 
after the edge saturation 56 , depending on the radicals or groups. 
Accordingly, the edge physicochemical amount a 2D is suppressed, 
leading to the change in E g (D) in an decreasing order with R = null, 
OH, F and H as observed, Thanks to it, the opening suppression 
originated from the edge saturation can be understood associated 
with the edge physicochemical nature. 

In light of equations 2.1 and 2.2, moreover, Figure 6 gives a com- 
parison of D necessary to achieve E g (D) = 0.1 eV (denoted as D 0 leV ) 
and E g (D) = 0.4 eV (marked with D 0 4eV ) between naked and satu- 
rated GNRs, which are specified by E c of edge-C atoms along the x- 
axis. From naked GNRs to saturated GNRs with R = OH, F and H, 
interestingly, D decreases as E c rises, while the D 0AeV curve goes 
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Figure 5 | \E g (D) as the function of D for disordered GQDs and GNRs 
between the cases with the saturated and naked edges, where \E g (D) = 
E g (D) R — E g (D) nu \[ with R = H, F and OH. The corresponding values of 

a 2D> a 2D an d E>2 D for calculation are referenced to those adopted in 
Figures 2 and 3. 



Figure 6 | A plot of D necessary to achieve E g (D) = 0.1 eV and 
E g (D) = 0.4 eV for disordered GNRs with or without the edge saturation 
specified by E c of edge-C atoms. The values of lQ,E*, a 0 ^, a^, and 
necessary for the calculation are referenced to those adopted in Figures 2 
and 3. 
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Figure 7 | E g (D) as the function of x R at D = 1 .6 nm for disordered GQDs 
and GNRs with (A) R = H, (B) R = F and (C) R = OH using equations 2.1 
and 2.2, where the corresponding values of E° c > a^, a^ D and D° 2T) are 
referenced to those adopted in Figures 2 and 3. 

below the D 0 leV curve. Relative to the naked case, small D 0 leV and 
A).4ev are strictly required for GNR-R, decreasing in the order of 
R = OH, F and H. Owing to these, one can see that, if a certain E g (D) 
value is devised for a graphene flake, the necessary D will be affected 
by the edge physicochemical nature associated with the saturation. 

In published papers, special interest has been attracted how to tune 
the gap of graphene. One important way for it is to control the edge 
states or nature of graphene flakes via the edge saturation with 
different chemical radicals or groups 31 . In ZZ-GNRs, the edge 
modifications can break the spin degeneracy. This will lead to a 
semiconducting-metal transition or a half- semiconducting state with 
the two spin channels having different band gaps, or a spin-polarized 
half- semiconducting state with the spins in the valence and conduc- 
tion bands being oppositely polarized. About AC- GNRs, the edge 
saturation gives electronic states a few eV away from the Fermi level 
and does not significantly affect their bandgap. Based on it 31 , hence, it 
is reasonable that E g (D) of disordered GQDs or GNRs can be modu- 
lated quantitatively by the edge saturation. In view of our predictions, 
since there exists an obvious discrepancy of E g (D) between naked and 
saturated GQDs or GNRs, such a modulation can be realized by 
controlling the edge coverage upon the saturation, which can be 
performed through eliminating only a certain percentage of dangling 
bonds in a uniform way. Given that there exists a linear relationship 
between oc 2 d and the edge coverage, one has oc 2 d = a 2D^ ~ x r) + 
^2D"^R f° r partially saturated GQDs and GNRs, where x R denotes the 
percentage of saturated edge-C atoms from 0 (the full edge-naked 
state) to 1 (the complete edge-saturated state). Inserting the above 
oc 2 d expression into equation 2.2, E g (D) as the function of x R at D = 
1.6 nm was plotted in Figure 7 using equation 2.1 with R = H in (A), 
R = F in (B) and R = OH in (C). E g (D) decreases continuously as x R 
rises. On increasing x R from 0 to 1, E g (D) declines from 1.70 eV to 
0.93 eV for GQDs and from 0.84 eV to 0.41 eV for GNRs in (A) with 
R = H, from 1.70 eV to 1.12 eV for GQDs and from 0.84 eV to 0.50 



eV for GNRs in (B) with R = F, and from 1.70 eV to 1.18 eV for 
GQDs and from 0.84 eV to 0.53 eV for GNRs in (C) with R = OH. 
This suggests that E g (D) of disordered GQDs and GNRs can be 
adjusted subtly by changing x R , beneficial to the design and engin- 
eering of graphene-based nanoelectronic devices. 

Methods 

General consideration on the BOs in disordered graphene flakes. Although the 
dimensionality of 2D graphene differs from 3D semiconductors, there exists CI for 
edge-C atoms of graphene. By reason of such a similarity, the coherent energy of 
edge-C atoms will be decreased definitely. It is natural that the nearly- free-electron 
approach can also be adopted in predicting the BOs in graphene. In view of this, the 
decrease in E c of edge-C atoms will play an essential role in influencing the crystalline 
field of graphene, and the equivalence between the A and B sublattices will resultantly 
be broken, leading to the BOs in graphene. 

According to the tight-binding approach 57 , electrons are viewed as occupying the 
standard orbitals of their constituent atoms, and then hopping between atoms during 
conduction. Mathematically, this is represented as a hopping or transfer integral 
between neighboring atoms, which can be viewed as the physical principle that creates 
electron bands in crystalline materials, due to overlapping between atomic orbitals. In 
virtue of it, the BOs in graphene induced by the upward shift of conduction band and 
downward shift of valence band of the n bonds are associated to the hopping energy 
E h between n orbitals with E h = 2.96 eV 57 . In light of these, with reference to 3D 
nano semiconductors 44 ' 45 , E g (D) of graphene flakes is given as, 

£ g (D)=[l-E c (D) 2D /£ c (cx))]£ h . (2.1) 

To obtain E g (D) with equation 2.1, £ c (D) 2D AE c (o°) function should be explored for 
2D flakes, which can be referenced with the 3D case 44 ' 45 ' 49,50 , where free surfaces of 3D 
nanocrystals with broken bonds are substituted by free edges of 2D flakes. One thus 
has, 

E C (D) 2D /E c (cv) = exp [- (a 2D - 1)/(D/D° D - l)] . (2.2) 

In equation 2.2, oc 2D is an edge physicochemical amount, while D 0 m denotes the 
critical size of 2D materials. To explore £ C (D) 2D , these two amounts should be 
developed by considering the roles played by the dimension, geometric shape and/or 
the edge saturation. Referenced by the 3D case 44 ' 45 ' 49,50 , the physicochemical amount 
oc 2D in equation 2.2 for 2D graphene flakes is decided by the edge nature relative to the 
interior case, and one thus has, 

0t2D = CTedge(^)7 CT in-2D(D) 2 , (2.3) 

where the subscripts edge and in-2D mean the edge and interior atoms of 2D graphene 
flakes. In principle, a 2D of naked graphene flakes marked with a 2D is decided by the CI 
of C edge-atoms. Since the CI of edge-C atoms is eliminated upon the edge saturation, 
in contrast, a different a 2D will be resulted from it for saturated graphene flakes 
denoted with o^ D . On the other hand, D Q 2T) is one amount related to the dimensionality 
or the geometrical shape of graphene flakes, which can be principally given as, 



D/D° 2D = s/lh, 



(2.4) 



where s/l is the area/edge ratio. Associated to it, the D- dependences of E g (D) should be 
relevant to the l/s ratio of graphene flakes. 

DFT approach to electronic structure of DO-GNRs. The electronic structure of 
DO-GNRs was investigated with the DFT approach. DO-GNRs are 1.2 nm wide or so 
with the AC or ZZ skeleton along the periodic axis, referred to as DO-GNRs(AC) or 
DO-GNRs(ZZ). Moreover, the two edges of both DO-GNRs(AC) and DO- 
GNRs(ZZ) are composed of the mixed AC or ZZ structures and kept symmetric. For 
DO-GNRs(AC), the supercell includes 12 zigzag C=C bonds in the periodic 
direction, as shown in Figure 4(a)-(c), where the edges are saturated by H, F and OH, 
respectively. For DO-GNRs(ZZ), the supercell has 12 armchair C=C bonds in the 
periodic direction [see Figure 4(d)], where the edge is saturated by H atoms. 

Spin-polarized DFT calculations are performed using DMOL3 code with GGA- 
PBE as the exchange correlation function. All electron core treatment and double 
numeric plus polarization (DNP) basis set are adopted. The Brillouin zone is sampled 
by 5 X 1 X 1 k-points and the real space global cutoff radius is set to be 5.0 A for the 
geometry and electronic calculations. In addition, a uniform vacuum of 15 A per- 
pendicular to the GNR is employed to ensure no interactions between neighboring 
cells. These bring out the convergence tolerance of energy of 1 X 10" 5 Ha, maximum 
force of 0.002 Ha A" 1 and maximum displacement of 0.005 A, respectively. 
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